Load Data

library(gdata)
dataset <- read.delim("raw_data/Figure3D.txt", stringsAsFactors = FALSE)

dataset$genotype <- factor(gsub(" ","",dataset$genotype))
dataset$experiment <- factor(rep(paste0("exp", 1:(length(dataset$genotype)/length(levels(dataset$genotype)))), each= length(levels(dataset$genotype))))

dataset$UID <- factor(paste(dataset$experiment, dataset$genotype))

kable(dataset, row.names = F)
genotype NT X2h X4h X6h experiment UID
WT 8.366667 35.60000 29.56667 21.23333 exp1 exp1 WT
PARP1KO 6.700000 57.23333 46.26667 45.96667 exp1 exp1 PARP1KO
ALC1KO 6.533333 56.96667 59.76667 55.66667 exp1 exp1 ALC1KO
ALC1KOPARP1KO 6.400000 55.26667 43.00000 46.40000 exp1 exp1 ALC1KOPARP1KO
WT 7.300000 36.13333 29.73333 23.00000 exp2 exp2 WT
PARP1KO 7.966667 48.96667 49.60000 43.93333 exp2 exp2 PARP1KO
ALC1KO 10.600000 54.16667 59.83333 57.93333 exp2 exp2 ALC1KO
ALC1KOPARP1KO 8.300000 53.43333 46.23333 42.43333 exp2 exp2 ALC1KOPARP1KO
WT 7.233333 34.06667 35.43333 20.66667 exp3 exp3 WT
PARP1KO 7.233333 39.76667 41.26667 33.56667 exp3 exp3 PARP1KO
ALC1KO 7.233333 54.46667 52.10000 55.13333 exp3 exp3 ALC1KO
ALC1KOPARP1KO 6.433333 42.43333 42.90000 35.96667 exp3 exp3 ALC1KOPARP1KO
library(reshape2)
# reshape to long format
dataset <- melt(dataset, variable.name = "treatment", value.name = "Outcome")

dataset$genotype <- relevel(dataset$genotype, ref = "WT")
dataset$UID <- relevel(dataset$UID, ref = "exp1 WT")

dataset$time <- gsub("NT","0",dataset$treatment)
dataset$time <- gsub("X|h","",dataset$time)
dataset$time <- as.integer(dataset$time)



kable(dataset, row.names = F)
genotype experiment UID treatment Outcome time
WT exp1 exp1 WT NT 8.366667 0
PARP1KO exp1 exp1 PARP1KO NT 6.700000 0
ALC1KO exp1 exp1 ALC1KO NT 6.533333 0
ALC1KOPARP1KO exp1 exp1 ALC1KOPARP1KO NT 6.400000 0
WT exp2 exp2 WT NT 7.300000 0
PARP1KO exp2 exp2 PARP1KO NT 7.966667 0
ALC1KO exp2 exp2 ALC1KO NT 10.600000 0
ALC1KOPARP1KO exp2 exp2 ALC1KOPARP1KO NT 8.300000 0
WT exp3 exp3 WT NT 7.233333 0
PARP1KO exp3 exp3 PARP1KO NT 7.233333 0
ALC1KO exp3 exp3 ALC1KO NT 7.233333 0
ALC1KOPARP1KO exp3 exp3 ALC1KOPARP1KO NT 6.433333 0
WT exp1 exp1 WT X2h 35.600000 2
PARP1KO exp1 exp1 PARP1KO X2h 57.233333 2
ALC1KO exp1 exp1 ALC1KO X2h 56.966667 2
ALC1KOPARP1KO exp1 exp1 ALC1KOPARP1KO X2h 55.266667 2
WT exp2 exp2 WT X2h 36.133333 2
PARP1KO exp2 exp2 PARP1KO X2h 48.966667 2
ALC1KO exp2 exp2 ALC1KO X2h 54.166667 2
ALC1KOPARP1KO exp2 exp2 ALC1KOPARP1KO X2h 53.433333 2
WT exp3 exp3 WT X2h 34.066667 2
PARP1KO exp3 exp3 PARP1KO X2h 39.766667 2
ALC1KO exp3 exp3 ALC1KO X2h 54.466667 2
ALC1KOPARP1KO exp3 exp3 ALC1KOPARP1KO X2h 42.433333 2
WT exp1 exp1 WT X4h 29.566667 4
PARP1KO exp1 exp1 PARP1KO X4h 46.266667 4
ALC1KO exp1 exp1 ALC1KO X4h 59.766667 4
ALC1KOPARP1KO exp1 exp1 ALC1KOPARP1KO X4h 43.000000 4
WT exp2 exp2 WT X4h 29.733333 4
PARP1KO exp2 exp2 PARP1KO X4h 49.600000 4
ALC1KO exp2 exp2 ALC1KO X4h 59.833333 4
ALC1KOPARP1KO exp2 exp2 ALC1KOPARP1KO X4h 46.233333 4
WT exp3 exp3 WT X4h 35.433333 4
PARP1KO exp3 exp3 PARP1KO X4h 41.266667 4
ALC1KO exp3 exp3 ALC1KO X4h 52.100000 4
ALC1KOPARP1KO exp3 exp3 ALC1KOPARP1KO X4h 42.900000 4
WT exp1 exp1 WT X6h 21.233333 6
PARP1KO exp1 exp1 PARP1KO X6h 45.966667 6
ALC1KO exp1 exp1 ALC1KO X6h 55.666667 6
ALC1KOPARP1KO exp1 exp1 ALC1KOPARP1KO X6h 46.400000 6
WT exp2 exp2 WT X6h 23.000000 6
PARP1KO exp2 exp2 PARP1KO X6h 43.933333 6
ALC1KO exp2 exp2 ALC1KO X6h 57.933333 6
ALC1KOPARP1KO exp2 exp2 ALC1KOPARP1KO X6h 42.433333 6
WT exp3 exp3 WT X6h 20.666667 6
PARP1KO exp3 exp3 PARP1KO X6h 33.566667 6
ALC1KO exp3 exp3 ALC1KO X6h 55.133333 6
ALC1KOPARP1KO exp3 exp3 ALC1KOPARP1KO X6h 35.966667 6

Plot Data

library(ggplot2)

# quadratic
ggplot(dataset, aes(x=time, y=Outcome)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE, aes(colour=genotype)) +
        geom_point(aes(colour=genotype, shape=experiment), size=2) +        
        facet_grid(. ~ genotype) +
        xlab(label = "Time post 2 Gy irradiation") +
        scale_shape_manual(values=15:20) +
        scale_color_manual(values=c("#000000","#000080","#808080","#800000"))

# cubic
ggplot(dataset, aes(x=time, y=Outcome)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE, aes(colour=genotype)) +
        geom_point(aes(colour=genotype, shape=experiment), size=2) +        
        facet_grid(. ~ genotype) +
        xlab(label = "Time post 2 Gy irradiation") +
        scale_shape_manual(values=15:20) +
        scale_color_manual(values=c("#000000","#000080","#808080","#800000"))

dataplot <- aggregate(Outcome ~ genotype+treatment, data = dataset, mean)
dataplot$se <- aggregate(Outcome ~ genotype+treatment, data = dataset, FUN = function(x){ sd(x)/sqrt(length(x))})$Outcome
dataplot$treatment <- factor(gsub("X|.nM","", dataplot$treatment)) 

dataplot$genotype <- factor(dataplot$genotype, levels = c("WT","PARP1KO","ALC1KO","ALC1KOPARP1KO"))
dataplot$treatment <- relevel(dataplot$treatment, ref = "NT")


ggplot(dataplot, aes(x=treatment, y=Outcome, fill = genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_bar(stat="identity", position=position_dodge()) +
        geom_errorbar(aes(ymin=Outcome-se, ymax=Outcome+se),width=.2, position=position_dodge(.9)) +
        ylab(label = paste0(expression("\u03B3"),"H2AX foci per cell")) +
        xlab(label = "Time post 2 Gy irradiation") +
        ylim(0, 70) +
        scale_fill_manual(values=c("#000000","#800000","#000080","#808080"))

library(Cairo)

cairo_pdf("Figure3D_v1.pdf", width = 7, height = 4, family = "Arial")

ggplot(dataplot, aes(x=treatment, y=Outcome, fill = genotype)) + 
        theme_bw() +
        theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(), 
              axis.line = element_line(colour = "black"), text = element_text(size=14),
              panel.border = element_blank(), panel.background = element_blank()) +
        geom_bar(stat="identity", position=position_dodge()) +
        geom_errorbar(aes(ymin=Outcome-se, ymax=Outcome+se),width=.2, position=position_dodge(.9)) +
        ylab(label = paste0(expression("\u03B3"),"H2AX foci per cell")) +
        xlab(label = "Time post 2 Gy irradiation (h)") +
        ylim(0, 70) +
        scale_fill_manual(values=c("#000000","#800000","#000080","#808080"))

dev.off()
## quartz_off_screen 
##                 2
dataset$genotype <- factor(dataset$genotype, levels = c("WT","PARP1KO","ALC1KO","ALC1KOPARP1KO"))
dataset$treatment <- factor(gsub("X","", dataset$treatment)) 
dataset$treatment <- relevel(dataset$treatment, ref = "NT")



cairo_pdf("Figure3D_v2.pdf", width = 7, height = 4, family = "Arial")

ggplot(dataplot, aes(x=treatment, y=Outcome, fill = genotype)) + 
        theme_bw() +
        theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(), 
              axis.line = element_line(colour = "black"), text = element_text(size=14),
              panel.border = element_blank(), panel.background = element_blank()) +
        geom_bar(stat="identity", position=position_dodge(), colour="black") +
        ylab(label = paste0(expression("\u03B3"),"H2AX foci per cell")) +
        xlab(label = "Time post 2 Gy irradiation (h)") +
        ylim(0, 70) +
        scale_fill_manual(values=c("#00000080","#80000080","#00008080","#80808080")) +
        geom_jitter(data = dataset, cex=1, position=position_dodge(0.9),  aes(colour=genotype)) +
        scale_color_manual(values=c("#000000","#800000","#000080","#808080")) +
        geom_errorbar(aes(ymin=Outcome-se, ymax=Outcome+se),width=.2, position=position_dodge(.9)) 


cairo_pdf("Figure3D_v3.pdf", width = 7, height = 4, family = "Arial")

ggplot(dataset, aes(x=treatment, y=Outcome, color=genotype)) +
        theme_bw() +
        theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(), 
              axis.line = element_line(colour = "black"), text = element_text(size=14),
              panel.border = element_blank(), panel.background = element_blank()) +
        geom_jitter(cex=1, position=position_dodge(1)) +
        stat_summary(fun.data=mean_se, fun.args = list(mult=1), geom="errorbar", width=0.25, aes(colour=genotype), position=position_dodge(1)) +
        stat_summary(fun.y=mean, geom="crossbar", width=0.5, aes(colour=genotype), position=position_dodge(1)) +
        ylab(label = paste0(expression("\u03B3"),"H2AX foci per cell")) +
        xlab(label = "Time post 2 Gy irradiation (h)") +
        ylim(0, 70) +
        scale_fill_manual(values=c("#00000080","#80000080","#00008080","#80808080")) +
        scale_color_manual(values=c("#000000","#800000","#000080","#808080"))


dev.off()
## quartz_off_screen 
##                 2

Linear Model

library(MASS)
library(DHARMa)
library(lme4)
library(lmerTest)
library(bbmle)
fit1 <- lm(Outcome ~ experiment+time*genotype, data = dataset)
print(summary(fit1))
## 
## Call:
## lm(formula = Outcome ~ experiment + time * genotype, data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -17.090 -11.403  -2.183  11.781  25.076 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                19.66444    7.29796   2.695   0.0104 *
## experimentexp2             -0.08542    4.92028  -0.017   0.9862  
## experimentexp3             -4.06458    4.92028  -0.826   0.4139  
## time                        1.91556    1.79663   1.066   0.2931  
## genotypePARP1KO             2.63111    9.50689   0.277   0.7835  
## genotypeALC1KO              3.95889    9.50689   0.416   0.6794  
## genotypeALC1KOPARP1KO       2.88556    9.50689   0.304   0.7631  
## time:genotypePARP1KO        3.01556    2.54082   1.187   0.2427  
## time:genotypeALC1KO         5.40444    2.54082   2.127   0.0400 *
## time:genotypeALC1KOPARP1KO  2.95111    2.54082   1.161   0.2527  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 13.92 on 38 degrees of freedom
## Multiple R-squared:  0.5491, Adjusted R-squared:  0.4423 
## F-statistic: 5.141 on 9 and 38 DF,  p-value: 0.0001444
cat("AIC: ", AIC(fit1))
## AIC:  399.7809
simres <- simulateResiduals(fittedModel = fit1)
plot(simres)

fit2 <- lm(Outcome ~ experiment+poly(time,2)*genotype, data = dataset)
print(summary(fit2))
## 
## Call:
## lm(formula = Outcome ~ experiment + poly(time, 2) * genotype, 
##     data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.4611  -4.2647  -0.0672   3.4001  13.5978 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           25.41111    2.28095  11.141 6.84e-13 ***
## experimentexp2                        -0.08542    2.28095  -0.037   0.9703    
## experimentexp3                        -4.06458    2.28095  -1.782   0.0837 .  
## poly(time, 2)1                        29.67566   12.90300   2.300   0.0277 *  
## poly(time, 2)2                       -65.08662   12.90300  -5.044 1.51e-05 ***
## genotypePARP1KO                       11.67778    2.63381   4.434 9.19e-05 ***
## genotypeALC1KO                        20.17222    2.63381   7.659 6.66e-09 ***
## genotypeALC1KOPARP1KO                 11.73889    2.63381   4.457 8.59e-05 ***
## poly(time, 2)1:genotypePARP1KO        46.71679   18.24760   2.560   0.0151 *  
## poly(time, 2)2:genotypePARP1KO       -14.43376   18.24760  -0.791   0.4344    
## poly(time, 2)1:genotypeALC1KO         83.72529   18.24760   4.588 5.83e-05 ***
## poly(time, 2)2:genotypeALC1KO        -18.16729   18.24760  -0.996   0.3265    
## poly(time, 2)1:genotypeALC1KOPARP1KO  45.71842   18.24760   2.505   0.0172 *  
## poly(time, 2)2:genotypeALC1KOPARP1KO -14.20282   18.24760  -0.778   0.4418    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.452 on 34 degrees of freedom
## Multiple R-squared:  0.9133, Adjusted R-squared:  0.8801 
## F-statistic: 27.55 on 13 and 34 DF,  p-value: 2.899e-14
cat("AIC: ", AIC(fit2))
## AIC:  328.6398
simres <- simulateResiduals(fittedModel = fit2)
plot(simres)

fit3 <- lm(Outcome ~ experiment+poly(time,3)*genotype, data = dataset)
print(summary(fit3))
## 
## Call:
## lm(formula = Outcome ~ experiment + poly(time, 3) * genotype, 
##     data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -6.2076 -1.9667  0.0132  1.7253  7.1944 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           25.41111    1.21180  20.970  < 2e-16 ***
## experimentexp2                        -0.08542    1.21180  -0.070 0.944274    
## experimentexp3                        -4.06458    1.21180  -3.354 0.002170 ** 
## poly(time, 3)1                        29.67566    6.85499   4.329 0.000154 ***
## poly(time, 3)2                       -65.08662    6.85499  -9.495 1.51e-10 ***
## poly(time, 3)3                        19.41656    6.85499   2.832 0.008176 ** 
## genotypePARP1KO                       11.67778    1.39927   8.346 2.58e-09 ***
## genotypeALC1KO                        20.17222    1.39927  14.416 5.01e-15 ***
## genotypeALC1KOPARP1KO                 11.73889    1.39927   8.389 2.31e-09 ***
## poly(time, 3)1:genotypePARP1KO        46.71679    9.69443   4.819 3.88e-05 ***
## poly(time, 3)2:genotypePARP1KO       -14.43376    9.69443  -1.489 0.146959    
## poly(time, 3)3:genotypePARP1KO        13.65011    9.69443   1.408 0.169399    
## poly(time, 3)1:genotypeALC1KO         83.72529    9.69443   8.636 1.24e-09 ***
## poly(time, 3)2:genotypeALC1KO        -18.16729    9.69443  -1.874 0.070701 .  
## poly(time, 3)3:genotypeALC1KO         13.13372    9.69443   1.355 0.185606    
## poly(time, 3)1:genotypeALC1KOPARP1KO  45.71842    9.69443   4.716 5.19e-05 ***
## poly(time, 3)2:genotypeALC1KOPARP1KO -14.20282    9.69443  -1.465 0.153310    
## poly(time, 3)3:genotypeALC1KOPARP1KO  22.06740    9.69443   2.276 0.030131 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.427 on 30 degrees of freedom
## Multiple R-squared:  0.9784, Adjusted R-squared:  0.9662 
## F-statistic: 79.96 on 17 and 30 DF,  p-value: < 2.2e-16
cat("AIC: ", AIC(fit3))
## AIC:  269.9136
simres <- simulateResiduals(fittedModel = fit3)
plot(simres)

Compare Results

ICtab(fit1,fit2,fit3,
      base=T)
##      AIC   dAIC  df
## fit3 269.9   0.0 19
## fit2 328.6  58.7 15
## fit1 399.8 129.9 11

Final Result

fit <- fit3

output <- coef(summary(fit))
output <- output[grep("time", rownames(output)),]

rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype",  paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1],  sep = " in " )


# suggested result table
kable(output, row.names = T)
Estimate Std. Error t value Pr(>|t|)
time1 in WT 29.67566 6.854994 4.329057 0.0001535
time2 in WT -65.08662 6.854994 -9.494774 0.0000000
time3 in WT 19.41656 6.854994 2.832469 0.0081758
time1: WT vs. PARP1KO 46.71679 9.694426 4.818933 0.0000388
time2: WT vs. PARP1KO -14.43376 9.694426 -1.488872 0.1469593
time3: WT vs. PARP1KO 13.65011 9.694426 1.408038 0.1693989
time1: WT vs. ALC1KO 83.72529 9.694426 8.636437 0.0000000
time2: WT vs. ALC1KO -18.16729 9.694426 -1.873993 0.0707005
time3: WT vs. ALC1KO 13.13372 9.694426 1.354770 0.1856062
time1: WT vs. ALC1KOPARP1KO 45.71842 9.694426 4.715949 0.0000519
time2: WT vs. ALC1KOPARP1KO -14.20282 9.694426 -1.465050 0.1533098
time3: WT vs. ALC1KOPARP1KO 22.06740 9.694426 2.276298 0.0301313
write.table(output, file = "Figure3D_Stats_Ref_WT.txt", quote = F, sep = "\t", row.names = T, col.names = NA)
# re-fit with ALC1KO reference
dataset$genotype <- relevel(dataset$genotype, ref = "ALC1KO")

fit <- lm(Outcome ~ experiment+poly(time,3)*genotype, data = dataset)

output <- coef(summary(fit))
output <- output[grep("time", rownames(output)),]

rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype",  paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1],  sep = " in " )


# suggested result table
kable(output, row.names = T)
Estimate Std. Error t value Pr(>|t|)
time1 in ALC1KO 113.4009524 6.854994 16.5428227 0.0000000
time2 in ALC1KO -83.2539088 6.854994 -12.1450008 0.0000000
time3 in ALC1KO 32.5502734 6.854994 4.7484028 0.0000474
time1: ALC1KO vs. WT -83.7252933 9.694426 -8.6364367 0.0000000
time2: ALC1KO vs. WT 18.1672885 9.694426 1.8739933 0.0707005
time3: ALC1KO vs. WT -13.1337169 9.694426 -1.3547700 0.1856062
time1: ALC1KO vs. PARP1KO -37.0085075 9.694426 -3.8175039 0.0006289
time2: ALC1KO vs. PARP1KO 3.7335317 9.694426 0.3851215 0.7028645
time3: ALC1KO vs. PARP1KO 0.5163978 9.694426 0.0532675 0.9578719
time1: ALC1KO vs. ALC1KOPARP1KO -38.0068766 9.694426 -3.9204877 0.0004749
time2: ALC1KO vs. ALC1KOPARP1KO 3.9644719 9.694426 0.4089434 0.6854860
time3: ALC1KO vs. ALC1KOPARP1KO 8.9336816 9.694426 0.9215277 0.3641269
write.table(output, file = "Figure3D_Stats_Ref_ALC1KO.txt", quote = F, sep = "\t", row.names = T, col.names = NA)

ANOVA

fit3a <- lm(Outcome ~ experiment+poly(time,3)*genotype, data = dataset)
fit3b <- lm(Outcome ~ experiment+poly(time,3)+genotype, data = dataset)

anova(fit3a, fit3b)
## Analysis of Variance Table
## 
## Model 1: Outcome ~ experiment + poly(time, 3) * genotype
## Model 2: Outcome ~ experiment + poly(time, 3) + genotype
##   Res.Df     RSS Df Sum of Sq      F    Pr(>F)    
## 1     30  352.43                                  
## 2     39 1343.95 -9   -991.51 9.3778 1.297e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

New analysis

dataset$genotype <- relevel(dataset$genotype, ref = "WT")


for(i in 2:length(levels(dataset$treatment))){
        
        dataset_sub <- dataset[dataset$treatment == levels(dataset$treatment)[i],]
        
        for(g in seq_along(levels(dataset_sub$genotype))){
                
                
                cat("Time: ", levels(dataset$treatment)[i], "\n",
                    "Genotype: ", levels(dataset$genotype)[g], "\n")
                
                
                dataset_sub$genotype <- relevel(dataset_sub$genotype, ref = levels(dataset_sub$genotype)[g])

                fit_sub <- lm(Outcome ~ genotype, data = dataset_sub)
                
                print(summary(fit_sub))
                cat("AIC: ", AIC(fit_sub))
                
                
                simres <- simulateResiduals(fittedModel = fit_sub)
                plot(simres)
                
                if(i == 2 & g == 1){
                        output <- coef(summary(fit_sub))
                        output <- output[grep("genotype", rownames(output)),]
                        rownames(output) <- paste0(gsub("genotype",  paste0(" ",levels(dataset$genotype)[g], " vs. "), rownames(output)),
                                                   " in Time ", levels(dataset$treatment)[i])
                } else {
                        outtmp <- coef(summary(fit_sub))
                        outtmp <- outtmp[grep("genotype", rownames(outtmp)),]
                        rownames(outtmp) <- paste0(gsub("genotype",  paste0(" ",levels(dataset$genotype)[g], " vs. "), rownames(outtmp)),
                                                   " in Time ", levels(dataset$treatment)[i])
                        
                        output <- rbind(output,outtmp)
                        
                        if(g == length(levels(dataset_sub$genotype)) & i < length(levels(dataset$treatment))){
                                
                                output <- rbind(output, " ", colnames(output))
                        }
                }
                
        }
        
        
}
## Time:  2h 
##  Genotype:  WT 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8889 -1.0750  0.3222  2.0889  8.5778 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             35.267      3.266  10.797 4.78e-06 ***
## genotypeALC1KO          19.933      4.619   4.315  0.00256 ** 
## genotypePARP1KO         13.389      4.619   2.898  0.01994 *  
## genotypeALC1KOPARP1KO   15.111      4.619   3.271  0.01134 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.658 on 8 degrees of freedom
## Multiple R-squared:  0.7191, Adjusted R-squared:  0.6137 
## F-statistic: 6.826 on 3 and 8 DF,  p-value: 0.01349
## 
## AIC:  80.78097

## Time:  2h 
##  Genotype:  ALC1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8889 -1.0750  0.3222  2.0889  8.5778 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             55.200      3.266  16.899 1.52e-07 ***
## genotypeWT             -19.933      4.619  -4.315  0.00256 ** 
## genotypePARP1KO         -6.544      4.619  -1.417  0.19431    
## genotypeALC1KOPARP1KO   -4.822      4.619  -1.044  0.32705    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.658 on 8 degrees of freedom
## Multiple R-squared:  0.7191, Adjusted R-squared:  0.6137 
## F-statistic: 6.826 on 3 and 8 DF,  p-value: 0.01349
## 
## AIC:  80.78097

## Time:  2h 
##  Genotype:  PARP1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8889 -1.0750  0.3222  2.0889  8.5778 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             48.656      3.266  14.896 4.07e-07 ***
## genotypeALC1KO           6.544      4.619   1.417   0.1943    
## genotypeWT             -13.389      4.619  -2.898   0.0199 *  
## genotypeALC1KOPARP1KO    1.722      4.619   0.373   0.7190    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.658 on 8 degrees of freedom
## Multiple R-squared:  0.7191, Adjusted R-squared:  0.6137 
## F-statistic: 6.826 on 3 and 8 DF,  p-value: 0.01349
## 
## AIC:  80.78097

## Time:  2h 
##  Genotype:  ALC1KOPARP1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8889 -1.0750  0.3222  2.0889  8.5778 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       50.378      3.266  15.423 3.11e-07 ***
## genotypePARP1KO   -1.722      4.619  -0.373   0.7190    
## genotypeALC1KO     4.822      4.619   1.044   0.3270    
## genotypeWT       -15.111      4.619  -3.271   0.0113 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.658 on 8 degrees of freedom
## Multiple R-squared:  0.7191, Adjusted R-squared:  0.6137 
## F-statistic: 6.826 on 3 and 8 DF,  p-value: 0.01349
## 
## AIC:  80.78097

## Time:  4h 
##  Genotype:  WT 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.1333 -1.8861 -0.2444  2.5500  3.8889 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             31.578      2.084  15.154 3.56e-07 ***
## genotypeALC1KO          25.656      2.947   8.706 2.36e-05 ***
## genotypePARP1KO         14.133      2.947   4.796  0.00136 ** 
## genotypeALC1KOPARP1KO   12.467      2.947   4.230  0.00288 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.609 on 8 degrees of freedom
## Multiple R-squared:  0.9049, Adjusted R-squared:  0.8693 
## F-statistic: 25.39 on 3 and 8 DF,  p-value: 0.0001931
## 
## AIC:  69.99344

## Time:  4h 
##  Genotype:  ALC1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.1333 -1.8861 -0.2444  2.5500  3.8889 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             57.233      2.084  27.465 3.33e-09 ***
## genotypeWT             -25.656      2.947  -8.706 2.36e-05 ***
## genotypePARP1KO        -11.522      2.947  -3.910  0.00448 ** 
## genotypeALC1KOPARP1KO  -13.189      2.947  -4.475  0.00207 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.609 on 8 degrees of freedom
## Multiple R-squared:  0.9049, Adjusted R-squared:  0.8693 
## F-statistic: 25.39 on 3 and 8 DF,  p-value: 0.0001931
## 
## AIC:  69.99344

## Time:  4h 
##  Genotype:  PARP1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.1333 -1.8861 -0.2444  2.5500  3.8889 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             45.711      2.084  21.936 1.97e-08 ***
## genotypeALC1KO          11.522      2.947   3.910  0.00448 ** 
## genotypeWT             -14.133      2.947  -4.796  0.00136 ** 
## genotypeALC1KOPARP1KO   -1.667      2.947  -0.566  0.58721    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.609 on 8 degrees of freedom
## Multiple R-squared:  0.9049, Adjusted R-squared:  0.8693 
## F-statistic: 25.39 on 3 and 8 DF,  p-value: 0.0001931
## 
## AIC:  69.99344

## Time:  4h 
##  Genotype:  ALC1KOPARP1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.1333 -1.8861 -0.2444  2.5500  3.8889 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       44.044      2.084  21.136 2.64e-08 ***
## genotypePARP1KO    1.667      2.947   0.566  0.58721    
## genotypeALC1KO    13.189      2.947   4.475  0.00207 ** 
## genotypeWT       -12.467      2.947  -4.230  0.00288 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.609 on 8 degrees of freedom
## Multiple R-squared:  0.9049, Adjusted R-squared:  0.8693 
## F-statistic: 25.39 on 3 and 8 DF,  p-value: 0.0001931
## 
## AIC:  69.99344

## Time:  6h 
##  Genotype:  WT 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.5889 -1.0028  0.2167  1.9611  4.8111 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             21.633      2.511   8.616 2.55e-05 ***
## genotypeALC1KO          34.611      3.551   9.747 1.03e-05 ***
## genotypePARP1KO         19.522      3.551   5.498 0.000575 ***
## genotypeALC1KOPARP1KO   19.967      3.551   5.623 0.000497 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.349 on 8 degrees of freedom
## Multiple R-squared:  0.9231, Adjusted R-squared:  0.8942 
## F-statistic: 31.99 on 3 and 8 DF,  p-value: 8.358e-05
## 
## AIC:  74.46729

## Time:  6h 
##  Genotype:  ALC1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.5889 -1.0028  0.2167  1.9611  4.8111 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             56.244      2.511  22.400 1.67e-08 ***
## genotypeWT             -34.611      3.551  -9.747 1.03e-05 ***
## genotypePARP1KO        -15.089      3.551  -4.249  0.00280 ** 
## genotypeALC1KOPARP1KO  -14.644      3.551  -4.124  0.00333 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.349 on 8 degrees of freedom
## Multiple R-squared:  0.9231, Adjusted R-squared:  0.8942 
## F-statistic: 31.99 on 3 and 8 DF,  p-value: 8.358e-05
## 
## AIC:  74.46729

## Time:  6h 
##  Genotype:  PARP1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.5889 -1.0028  0.2167  1.9611  4.8111 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)            41.1556     2.5109  16.391 1.93e-07 ***
## genotypeALC1KO         15.0889     3.5509   4.249 0.002802 ** 
## genotypeWT            -19.5222     3.5509  -5.498 0.000575 ***
## genotypeALC1KOPARP1KO   0.4444     3.5509   0.125 0.903482    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.349 on 8 degrees of freedom
## Multiple R-squared:  0.9231, Adjusted R-squared:  0.8942 
## F-statistic: 31.99 on 3 and 8 DF,  p-value: 8.358e-05
## 
## AIC:  74.46729

## Time:  6h 
##  Genotype:  ALC1KOPARP1KO 
## 
## Call:
## lm(formula = Outcome ~ genotype, data = dataset_sub)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.5889 -1.0028  0.2167  1.9611  4.8111 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      41.6000     2.5109  16.568 1.78e-07 ***
## genotypePARP1KO  -0.4444     3.5509  -0.125 0.903482    
## genotypeALC1KO   14.6444     3.5509   4.124 0.003325 ** 
## genotypeWT      -19.9667     3.5509  -5.623 0.000497 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.349 on 8 degrees of freedom
## Multiple R-squared:  0.9231, Adjusted R-squared:  0.8942 
## F-statistic: 31.99 on 3 and 8 DF,  p-value: 8.358e-05
## 
## AIC:  74.46729

# suggested result table
kable(output, row.names = T)
Estimate Std. Error t value Pr(>|t|)
WT vs. ALC1KO in Time 2h 19.9333333366666 4.61941688316092 4.31511895133979 0.00256288937859708
WT vs. PARP1KO in Time 2h 13.38888889 4.61941688316092 2.89839372125219 0.0199407600005576
WT vs. ALC1KOPARP1KO in Time 2h 15.11111111 4.61941688316092 3.27121614961496 0.0113355916961444
ALC1KO vs. WT in Time 2h -19.9333333366666 4.61941688316091 -4.31511895133979 0.00256288937859707
ALC1KO vs. PARP1KO in Time 2h -6.54444444666667 4.61941688316091 -1.41672523008759 0.194305008686908
ALC1KO vs. ALC1KOPARP1KO in Time 2h -4.82222222666667 4.61941688316091 -1.04390280172483 0.32704543157344
PARP1KO vs. ALC1KO in Time 2h 6.54444444666666 4.61941688316092 1.41672523008759 0.194305008686909
PARP1KO vs. WT in Time 2h -13.38888889 4.61941688316092 -2.89839372125219 0.0199407600005575
PARP1KO vs. ALC1KOPARP1KO in Time 2h 1.72222222 4.61941688316092 0.372822428362763 0.718964761240989
ALC1KOPARP1KO vs. PARP1KO in Time 2h -1.72222222 4.61941688316091 -0.372822428362763 0.71896476124099
ALC1KOPARP1KO vs. ALC1KO in Time 2h 4.82222222666667 4.61941688316092 1.04390280172483 0.32704543157344
ALC1KOPARP1KO vs. WT in Time 2h -15.11111111 4.61941688316091 -3.27121614961496 0.0113355916961444
Estimate Std. Error t value Pr(>|t|)
WT vs. ALC1KO in Time 4h 25.6555555566666 2.94700098560667 8.70564878735025 2.36474253843872e-05
WT vs. PARP1KO in Time 4h 14.1333333366667 2.94700098560667 4.79583597212716 0.00136277121174079
WT vs. ALC1KOPARP1KO in Time 4h 12.4666666666667 2.94700098560667 4.23028927630314 0.00287522676347967
ALC1KO vs. WT in Time 4h -25.6555555566667 2.94700098560666 -8.70564878735025 2.3647425384387e-05
ALC1KO vs. PARP1KO in Time 4h -11.52222222 2.94700098560666 -3.90981281522309 0.0044822300200186
ALC1KO vs. ALC1KOPARP1KO in Time 4h -13.18888889 2.94700098560666 -4.47535951104711 0.00206847786110044
PARP1KO vs. ALC1KO in Time 4h 11.52222222 2.94700098560667 3.90981281522308 0.00448223002001862
PARP1KO vs. WT in Time 4h -14.1333333366667 2.94700098560667 -4.79583597212717 0.00136277121174079
PARP1KO vs. ALC1KOPARP1KO in Time 4h -1.66666667 2.94700098560667 -0.565546695824028 0.587208761244944
ALC1KOPARP1KO vs. PARP1KO in Time 4h 1.66666667 2.94700098560666 0.565546695824028 0.587208761244943
ALC1KOPARP1KO vs. ALC1KO in Time 4h 13.18888889 2.94700098560667 4.47535951104711 0.00206847786110044
ALC1KOPARP1KO vs. WT in Time 4h -12.4666666666667 2.94700098560666 -4.23028927630314 0.00287522676347967
Estimate Std. Error t value Pr(>|t|)
WT vs. ALC1KO in Time 6h 34.61111111 3.55089103646649 9.747162262811 1.02732215316808e-05
WT vs. PARP1KO in Time 6h 19.5222222233333 3.55089103646649 5.49783759142325 0.000575210034763358
WT vs. ALC1KOPARP1KO in Time 6h 19.9666666666667 3.55089103646649 5.62300179352606 0.000496757913056301
ALC1KO vs. WT in Time 6h -34.61111111 3.55089103646648 -9.74716226281101 1.02732215316807e-05
ALC1KO vs. PARP1KO in Time 6h -15.0888888866667 3.55089103646648 -4.24932467138775 0.00280172251246753
ALC1KO vs. ALC1KOPARP1KO in Time 6h -14.6444444433333 3.55089103646648 -4.12416046928495 0.00332508800054009
PARP1KO vs. ALC1KO in Time 6h 15.0888888866667 3.55089103646649 4.24932467138775 0.00280172251246755
PARP1KO vs. WT in Time 6h -19.5222222233333 3.55089103646649 -5.49783759142326 0.000575210034763353
PARP1KO vs. ALC1KOPARP1KO in Time 6h 0.44444444333333 3.55089103646649 0.125164202102805 0.903481679641302
ALC1KOPARP1KO vs. PARP1KO in Time 6h -0.444444443333329 3.55089103646648 -0.125164202102805 0.903481679641302
ALC1KOPARP1KO vs. ALC1KO in Time 6h 14.6444444433333 3.55089103646649 4.12416046928495 0.0033250880005401
ALC1KOPARP1KO vs. WT in Time 6h -19.9666666666667 3.55089103646648 -5.62300179352606 0.000496757913056298
write.table(output, file = "Figure3D_Stats_New_All.txt", quote = F, sep = "\t", row.names = T, col.names = NA)